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1. Introduction 

Let X G M^^"'^-^ be a tensor of order 3, sometimes named a three-way array or a 
three-mode data set. A rank 1 or a decomposed tensor is 

D = a(g)b(g)c, (1) 

where a G M^, b G R"^ and c G M'^, and ® is the tensor product, sometimes named 
also outer product. X can be expressed as a sum of decomposed tensors given in (1), 

r 

X = ^D„. (2) 

The rank of X is defined to be the minimal integer r, see for instance Kruskal (1977, 
1989). In data analysis, this implies that the rank of a three-way array is the smallest 
number of components that provide a perfect fit in Candecomp/Parafac (CP), see for 
instance, (Carroll and Chang, 1970, and Harshman, 1970). In statistics CP is con- 
sidered a natural extension of singular value decomposition or principal components 
analysis to three-way data. 

There is quite a literature concerning the value of maximal rank, generic rank 
or typical rank of three-way arrays in the area of statistics, algebraic complexity 
theory and algebraic geometry. Some references, among others, are: Ja' Ja' (1979), 
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Kruskal (1977, 1983, 1989), Strassen (1983), Ten Bcrgc (1991, 2000, 2004a, 2004b), 
Ten Berge and Kiers (1999), Ten Berge and Stegeman (2006), Comon and Ten Berge 
(2008), Biirgisser et al (1997), Catalisano et al (2002), Friedland (2008) and Abo et 
al. (2006). Friedland (2008) provides an up todate survey with some new results on 
the generic rank of three-way arrays. 
First, wc give the following 

Definition 1: A dataset is called generic if its elements are randomly generated 
from a continuous distribution. 

The generic and typical ranks are defined in the following way by Comon and Ten 
Berge (2008): Given that the rank of I x J x K arrays is bounded, one can partition 
the arrays according to the rank values. Generic rank is defined to be true almost 
everywhere; while typical ranks are associated with the rank values that occur with 
positive probability. So, if there is a single typical rank, then it may be called generic 
rank; that is, a generic rank is typical, but the converse is not true. 

Ten Berge (2000) classified three-way arrays into three classes: very tall, tall and 
compact. Let X G ^^^'^^^ be a tensor of order 3 with I > J > K. The array X is 
called very tall when / > KJ; X is tall when KJ — J < I < KJ — 1; X is compact 
when I < KJ — J. The generic rank of the very tall arrays is very well known 
and easiest to prove: it is KJ. Ten Berge (2000) showed that all tall three-way 
arrays have generic rank /; and the tallest among the compact arrays, that is when 
/ = KJ — J, have typical rank {/, / + 1} • Ten Berge and Stegeman (2006) provided 
some further results on the compact case. Friedland (2008) showed that: typical 
rank(12 x 4 x 4) >12, typical rank(ll x 4 x 4) >11, and typical rank(7 xJxK) > / for 
(/, J, K) — ((J— 1)^-1-1, J, J) when J > 2. These results are all based on mathematical 
proofs. However, the rank computation problem has also been approached from a 
numerical point of view: Gomon and ten Berge (2008) and Friedland (2008) applied 
Terracini's lemma, based on the numerical calculation of the maximal rank of the 
Jacobian matrix of (2), to obtain numerically the generic rank of some three-way 
arrays. The numerical method based on Terracini's lemma, when used to evaluate 
rank over R, gives the generic rank when the typical rank is single- valued, and the 
smallest typical rank value otherwise. 

Two well known facts are: a) There is no known method to calculate the rank of a 
given three-way dataset, Martin (2004, AIM tensor workshop); b) A three-way array 
over R may have a different rank than the same array considered over C, (Kruskal, 
1989). 

Wc shall be concerned by the numerical computation of the rank of a three-way 
array over R only. 

Computationally, the most primitive approach to the numerical evaluation of the 
typical rank of three-way arrays is based on the alternating least square (ALS) min- 
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imization algorithm: It is to run ALS many times to convergence on many generic 
three-way arrays of a given format, and to check whether or not the fit is perfect for 
a given number of components. But as a referee remarked, this approach has 2 prob- 
lems: First, we do not know how many three-way arrays of a given format to examine 
before a valid inference can be drawn. For instance, when 100000 arrays have been 
examined and all seem to have the same rank a, it does not follow that a is indeed 
the generic rank for that array format. After all, a different rank may occur with an 
extremely small yet positive probability. Second, the decision of when to terminate 
ALS is hazardous, because even if the residual sum of squares is, say, exp(— 32), this 
does not prove that it is zero; in fact, it may have zero as infimum. The present paper 
relieves us from both above mentioned problems: It offers a straightforward method 
of determining the rank of any given array over M, based on inspection of the number 
of real roots of a system of certain polynomial equations. 

The real solution set of a system of polynomial equations is called semi-algebraic 
set in real algebraic geometry, see Basu, Pollack and Roy (2006) or Friedland (2008). 
Semi-algebraic sets are open sets and are composed of a finite union of connected 
components, where each component is called a basic semi- algebraic set. The main 
problem can be reformulated as: For a given tensor X over R calculate the number 
of connected components where each component is characterized by a unique real 
rank value. Our numerical results will shed some light on this. The numerical results 
on simulated datasets will be obtained by computing the Grobncr bases using Maple 
12 of the system of polynomial equations characterizing the dataset. We note that 
generic datasets and random numbers are generated from integers between —99 and 
99. 

The paper is organized as follows. In section 2 we present the main lemma which 
provides a necessary and sufficient condition that a three-way array can be expressed 
as a sum of a fixed number of decomposed tensors. All results in this paper will 
be based on this lemma. In section 3, we show how the lemma can be applied to 
compute the rank of a generic tensor over R numerically for some cases. In section 4, 
we show another application of the lemma for the computation of rank for nongeneric 
particular datasets. In section 5, we show how the lemma can be applied to compute 
the rank of generic / symmetric J x J arrays, named INDSCAL arrays, over M. And 
finally in section 6 we conclude. 

2. Main Lemma 

Let X e M^x-^x-f^ be a three -way dataset and 2 < K < J < I. The lemma provides 
a necessary and sufficent condition that the tensor X can be expressed as a sum of / 
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decomposed tensors; that is 

I 

a=l 
I 

= 5Z a„ (8) b« c„, (3) 

where {aol a = 1,...,/} is a basis for R^, Cq, G M'^ and bo- G M"'. Note that if 
(3) is true, then rankpC) <I. We denote by G M.^^"^ the kth shce in X for 
k — 1, K.MVe note that (3) can be written as 

Xfe = ^Cfc^aa^bc, for k^l,...,K, 

a=l 

= AD(c,)B' for k^l,...,K, (4) 

where A = (ai aa-.a/) G R^^^ B = (bi ba-bj) G R-^^^ C = (c^) G R^^^ and 
D(c^) = Diag(cfc) G R^^^ is a diagonal matrix with diagonal elements Cka- Note that 

the vector G R''^ represents the A;th row of C. 
We consider the system of polynomial equations 

s^Xfc = Cfeab^ for A; = 1, and q; = 1, (5) 

where {Sa\ o; = 1, ...,/} is a basis for R^ and G R^, and b^ G R-^. We note that 
(5) can be written as 

S'Xfc = D(cjB' for A; = 1, K, (6) 

where S has columns s^. 

Lemma 1: (6) is a necessary and sufficient condition for (4). 
Proof. Let I = AS', then (4) is true if and only if (6) is true. 

Remark 1: a) To see if (5) is true, we solve the system of polynomial equations 

s'Xfc = Cfeb' for k^l,...,K, (7) 

for s G R^, b G R-^ and c G R^. 

b) Wc note that (7) has two indeterminacies: It can be rewritten as s'^l^k — Ck*^'* 
for k = 1, ...jK, , where for instance, s,. = As for any scalar A 7^ , Ck* = fJ^Ck for 
any scalar n ^ 0, and b* = \h/ fi. To eliminate these indeterminacies, hereafter, we 
fix 

ci — 1 and Si — 1. (8) 
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c) Theorem 2.4 of Friedland (2008) provides another characterization for (3) or 
(4): It states that each shce G span{a.i ^ bi, a/ (g) b/) and the rank(X) equals 
the minimal dimension of the span{a.i (g) bi, a/ (8) hi). 

d) The necessary condition, when rank(X) —I, which was shown to be also suf- 
ficient afterwards, was used many times by Ten Berge and his coworkers. Ten Berge 
(2000), Ten Berge (2004a), and Ten Berge, Sidiropoulos and Rocci (2004). 

3. Rank computation 
We shall suppose in the sequel that X G M^^-^^^ is a generic three-way array and 2 < 
K < J < I < KJ. Then we have the following well known inequality: rank(X) >/. 
We will check if X has rank /. By the Main Lemma , the tensor X has rank /, if for 
parameter vectors s G R^, b G M"^ and c G the system of polynomial equations 
(7) subject to (8) have / real solutions (cfea, b^, s^) for a — 1, such that the 
elements of the set {sq,| a — 1,...,/} is a basis for that is, (7) with (8) has / 
real isolated solutions. Let us see how can we know if this is true. The system of 
polynomial equations (7) with (8) is equivalent to 

s'(Xfc - CfeXi) = 0' for A; = 2, K. (9) 
So the number of equations, neg, in (9) is 

neq — {K — 1) J, 

and the number of degrees of freedom or the number of free variables , df, is 

df={I-l) + {K-l), (11) 

because of (8) there are {K — 1) free 's and (/ — 1) free Si 's. 

We are interested in the study of the number of solutions of (9) over R for generic 
data. We distinguish three cases named, minimal when neq = df, overdetermined 
when df < neq, and, underdetermined when df > neq. We note that Abo et al. 
(2006) also distinguished three cases that they named subabundant, superabundant 
and equiabundant: these were used for induction purposes. 

3.1. Case 1: Minimal System(neq = df). When I = {K - 1)(J - 1) + 1, 
neq = df, and the system (9) is called minimal. The number of real solutions is 
bounded; an upper bound is provided by Khovanskii's theorem, see Sturmfels (2002), 

Theorem 1 (Khovanskii): Consider n polynomials in n variables involving 
m distinct monomials in total. The number of isolated roots in the positive orthant 
(i?+)" of any such system is at most 2^^\n + 1)"*. 
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In our case n = neq = df = {K — 1)J and m = / — 1 + {K — — 1) = K{I — 1). 
The number of isolated roots in the positive orthant {R+y-^ of any such system is at 
most 2^^\df + 1)"* and m is the number of distinct monomials in the system (9). 
So (9) may or may not have / real isolated solutions. In case (9) has / real isolated 
solutions, then rank(X) = /; otherwise we embed it, which is discussed later on. 

Example 1: 1x1x2 arrays: neq = df = I 

This class of arrays is discussed in detail by Ten Berge (1991), who showed that 
the typical rank of such arrays is {/, / + 1}. To check if the rank of a generic 7x7x2 
array is 7, it suffices to solve (9), which reduces to finding the real roots of the 
determinantal equation det(K2 — C2X1) = 0. If det(K2 — C2X1) = has 7 real roots, 
then rank(X.) = 7, otherwise rank(%) =7 + 1. Simulation results for 5000 generic 
3x3x2 arrays produced one real root 51.76% and 3 real roots 48.24% of the time. 
So we can deduce that Pr(rank (3x3x2 array) = 3) ^ 48.24% and Pr(rank (3x3x2 
array) = 4) 51.76%. 

Example 2: 7 x J x 3 arrays with 7 = 2 J - 1: neq = rf/ = 2 J 
a) 5 x 3 x 3 arrays: neq = df = Q. This class of arrays is also discussed in Ten 
Berge (2004a), where Ten Berge showed that generic 5x3x3 arrays have either rank 
5 or rank 6 with positive probability. Further, he showed that a closed form solution 
for the case when the array has rank 5 corresponds to finding the number of real roots 
of a sixth degree polynomial equation: if there are 6 real roots, then the array has 
rank 5, otherwise its rank is 6. Table 1 displays the number of real roots obtained 
by solving the system (9) for 1000 simulated generic arrays. First, we note that 
the solution set of (9) always admitted 6 roots, as expected according to Ten Berge 
(2004a); further, the number of real solutions is an even number or zero. Second, 
Pr(rank (5x3x3 array) = 5)^ 6.8% and Pr(rank (5x3x3 array) = 6) 93.2%. 



Table 1: Simulation results for 1000 generic 5x3x3 arrays. 


real roots 





2 


4 


6 


counts 


47 


501 


384 


68 



b) 7 X 4 X 3 Eirrays: neq — df — 8. Table 2 displays the number of real roots 

obtained by solving the system (9) for 1000 simulated generic arrays. First, we note 
that the solution set of (9) always admitted 10 roots and the number of real solutions 
is an even number or zero. Second, Pr(rank (7x4x3 array) = 7) ~ 4.2%. 



Table 2: Simulation results for 1000 generic 7x4x3 arrays. 


real roots 





2 


4 


6 


8 


10 


counts 


16 


268 


456 


218 


40 


2 



c) 9 X 5 X 3 arrays: neq = df = 10. Table 3 displays the number of real roots 
obtained by solving the system (9) for 1000 simulated generic arrays. First, we note 
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that the solution set of (9) always admitted 15 roots and the number of real solutions 
is an odd number. Second, Pr(rank (9x5x3 array) = 9) ~ 6% and Pr(rank (9x5x3 
array) = 10) 94%. This latter result follows from Ten Berge (2000, Result 5) or 
see Example 5 describing tallest compact arrays, by embedding 9x5x3 arrays into 
10 X 5 X 3 arrays. 



Table 3: Simulation results for 1000 generic 9x5x3 arrays. 


real roots 


1 


3 


5 


7 


9 


11 


13 


counts 


34 


290 


404 


212 


51 


8 


1 



Example 3: / x J x 4 arrays with / = 3J — 2: neq = df = 3J 

Numerical computations showed that ^ {roots of 10 x 4 x 4 arrays) = 20; 7^ {roots 

of 13 X 5 X 4 arrays) = 35 and 7^ (roots of 16 x 6 x 4 arrays) = 56. Table 4 shows 

that Pr(rank (10 x 1 x 1 array) = 10) ^ 7.8%. 



Table 4: Simulation results for 1000 generic 10 x 4 x 4 Eirrays. 


real roots 





2 


4 


6 


8 


10 


12 


14 


counts 


2 


78 


284 


342 


216 


58 


14 


6 



Remark 2: a) To calculate a Grobner basis for (9) in Example 2 for Jx Jx3 arrays 
with / = 2 J — 1, we used pure lexicographic order given by the following sequence 
(si, C3, C2) of the free variables. In all cases the Grobner basis, denoted by Gp, 
consisted of {K — 1)J polynomials having the following form: 6*1(02) = 0, G2(c2, C3) = 
poly2(c2) + C3 = 0, G3^a{c2,Sa) = polyQ(c2) + Sa = foT a = 1,...,/ - 1. It is 
important to note that this particular form of the Grobner basis polynomials, Gfs, 
shows that the degree of the polynomial Gi(c2) = 0, denoted by dcgGi{c2), represents 
the number of roots of the system (9). An introduction to Grobner basis can be found 
in, among others, Gox et al. (2007). Example 6 show quite in detail the Grobner 
basis application to a generic array. 

b) The Maple 12 commands to do the computations in Example 2 are shown in 
Appendix 1. 

c) For 7x7x2 arrays and 7 > 2, (iet(X2 -C2X1) = Gi(c2) = 0, where Gi(c2) = 
is the first element of the Grobner basis. This phenomenon will be also seen for tallest 
compact arrays, see Examples 4, 5 and 6. 

A reviewer noted that the right hand side of (7) is a Segre variety, which is 
the image of the Segre map, S(^_i) The Segre map sends an clement of the 

projective space p(-^^^) x p('^^^) into P^-'~^. While the left hand side of (7) is a 
linear space of projective dimension 7 — 1 = {K — 1){J — 1). So, (7) , will represent the 
intersection of the linear space with the Segre map, and the number of intersections 
is the degree of the Segre variety given in (12), see for instance Harris (1992, p. 233). 
This result is summarized in the following 



Some Numerical Results on the Rank of Generic Three- Way Arrays over M 



8 



Theorem 2: Let / = (/^ - 1)( J - 1) + 1 and 2 < < J < J, then for generic 
data the number of roots (real or complex) of the polynomial system (9) is 

degG,{c2)=(^~^^_'[~^). (12) 



CorollEiry 1: For minimal systems and 3<K<J<I,I< degGi{c2). 
Proof: Let n — J — 1 and m — K — 1. We have to show that 



mn + 1 < 



(m + n)! 
~ n!m! 

It is true for m = 2. For m > 3. wc have 



for 2 < m < n. 



{m + n)\ 



n\m\ 



(n + l)(n + 2) (n + m-2) 



> 



m m — 1 3 
(n + m — 1) (n + m) 

2 i 



(n + m — 1) (n + m) 



So, it is sufficient to show that (n + m)(n + m — 1) > 2{mn + 1), which is easily seen 
to be true. 

Corollary 2: The typical rank of arrays with a minimal system have more than 
one rank value and the minimum attained value is /. 

Proof: The rank of a generic array with a minimal system is /, if the number of 
real roots of Gi(c2) is greater than or equal to /; otherwise its rank is greater than /. 

We note that Corollary 2 generalizes Friedland (2008), who showed that: typical 
rank(7 x J x X) > 7 for (7, J, K) = ((J - if + 1, J, J) when J > 2. 

3.2. Case 2: Underdetermined System(df > neq). When {K - 1)(J - 1) + 
2 < 7 < 7J, df > neq, and the system (9) is called underdetermined. The upper 
bound for the number of isolated roots of (9) is infinity; so (9) may or may not have 
7 real isolated solutions: So the attained minimum bound for the rank of a generic 
three-way array is, 6min = I- Before discussing two general classes studied in detail 
by Ten Berge (2000), we introduce some notation. 
The system (9) can be written as 



sT = S' [(X2 - C2X1), (X3 - C3X1), ...,Xk - C;^Xi] = 0', (13) 
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where the number of columns of the matrix T is 

ncoir = (K-1)J, (14) 
= neq; 

and the number of rows of F is 

nrowr=7, (15) 
and r is a matrix function of the parameters C2, c^. We also define 

nbil (16) 

which represents the minimal number of Ck parameters that can be specialized to 
make the system of polynomial equations (13) linear. In algebraic geometry, the 
replacement of variables by specific values is called specialization. 
Example 4: Tall arrays: df — neq > nbil 

These are arrays when {K — 1)J < I < KJ and I > J > K, whose generic rank 
is 7, as shown by Ten Berge (2000, Result 2). This impfies that (15) > (14), that 
is / > [K — 1)J, or, df — neq > nbil = K — 1, where nbil is given in (16). By 
assigning random values to the {K — l)c^.'s in (13), we reduce (13) to a system of 
linear equations, which will have a solution for any generic data; so (13) will admit / 
real and isolated solutions; from which we deduce that the generic rank of tall arrays 
is I. 

Example 5: Tallest compact arrays: ncolF =nrowr and > 3 

These are arrays when I = J{K — 1), I > J > K and K > ?>. Note that we 

exclude 1x1x2 arrays for / > 2 discussed in Example 1. Ten Berge (2000, Results 

3, 4 and 5) discussed this case. 

When I — {K — 1)J and X > 3, it implies that (14) = (15), that is, F is a square 

matrix. Solving (13) for 's for k — 2,...,K is equivalent to solving detiV) =0. 

K 

The leading monomial in det{r) =0 is ]^c^- If J is an odd integer, then (13) will 

k=2 

have infinite number of real solutions: Assign random continuous numbers to 's 
for /c = 3, i^, and solve for C2. This corresponds to Result 5 in Ten Berge (2000), 
which states: When / = J[K — 1) and I > J > K and K > 3 and J is odd, then 
the typical rank is I. If J is an even integer, then (13) may have infinite number 
of real solutions or finite number of real solutions or real solution: For instance 
for J = 4 and K = 3, the polynomial /(c2,C3) = Sc'^c'^ + 1 has real solution, the 
polynomial /(c2,C3) = ic'^c'^ — 1 has infinite number of real and distinct solutions, 
and the polynomial /(c2, C3) = 303(03 — 1) has a finite number of real solutions. Ten 
Berge (2000) specifically discussed the case of 8 x 4 x 3 arrays, where he stated that 
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typical rank of such arrays is {8, 9} and for randomly sampled data the rank of 9 
is extremely rare. Similarly, Friedland(2008, Th.7.2) showed that typical rank of 
12 X 4 X 4 arrays has more than one value. We conducted a limited simulation study 
on generic 8x4x3 and 12 x 4 x 4 arrays; and each time we got / real isolated 

solutions. The simulation study was done in the following way: For a generic dataset 
let /(c2, C3..., Ck) = det(r) =0; assign random values to the parameters C3..., Ck, then 
solve for C2. This shows that for generic data, when / = J{K — 1) and K > 3 the 
rank is / with very high probability. Also, see example 6. 

3.3. Example 6. We consider a simulated generic dataset of size 7x4x3 having 
the following three slices 

/ [-50,-38,-98,-93,-32,8,44] \ 

, [-22,-18,-77,-76,-74,69,92] 

1 [45, 87, 57, -72, -4, 99, -31] 

\ [-81,33,27,-2,27,29,67] / 
/ [99,-25,24,-61,31,25,50] \ 
, _ [60,51,65,-48,-50,94,10] 

2 [-95, 76, 86, 77, -80, 12, -16] 

\ [-20,-44,20,9,43,-2,-9] / 

/ [90,-82,29,52,42,-62,22] \ 
, [80,-70,70,-13,18,-33,14] 
^ [19, 41, -32, 82, -59, -68, 16] 

V [88,91,-1,72,12,-67,9] / 

Our aim is to find the rank of X, by representing it as in (6). This dataset has 
a minimal system of polynomial equations. We solve equation (9) via Grobner basis 
using the lexicograhic order {si, S2, S3, 84,85, sq, 03,02)- The first two polynomials of 
the Grobner basis are 

G'i(c2) = = 

-258797975083999058663603818114838724165583573114256294 

-1987946767932180125365724555441379125561037244553323732*C2 

-7447583055793225423658520296174635567495579387052082486*ci 

+18477292423934054741969006645285810935999768448664319668*ci 

+162868576676248184458245504688648649537661407605447407344*c^ 

+22324671325209561198922665813216562379249229294549244662*c^ 

-93044594774454916354246852601811640731920664188422515202*c^ 

+1034990365268175640254342731156079689724071145674294956746*4 

-1399215109838269848671482913176200716552825195700591390075*c| 

+155346700794650490501115016130585172314320583287574900147*c^ 
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+645072630378953757678717001000719217821315452777261680268*4° 
^2(02,03) = = 

23547338655791229204617338928506186026357940595072461145844743562575 
17224332611787380814907080082015469486345955534220731803730371551646 

24134092641462097641695580678924043455393253397537452244857196895154 

46461325491822959772758622053185481177195343628678109510431558799243 

57732629782484420615672656492158360322133109432630494581048689547 

-3097250590883846738286889878288351798304609282451856596410303546337 

36319437905844941144694584532179492363486873159383686238249006746894 

40160319769886619561748249254067684300497233927857246816557510045452 

29549030244795643979218209546091280718239589534262958680509337880753 

427578723122948642870383378005580770133313107270509303255691909320*C2 

-1903803609627374035621320978361604020106824484275906499213001947264 

95672033886734748658480731859038594412417938078895776898706276221196 

35127294306367325389388675757732402447319340359820738619607197620136 

52879058280068724250611547838816805478213287696041214145895488077794 

7039714261670007939310542108233323980869417122270525843479004885596*c2 

+8450512390839763624967161124974667624579807298666582130379782638299 

901903452364927022770162589839439882861673620370340791082608420840404 

311331236303835885178287071688967159857905276468573727498262705288023 

535530120799333316114432548548504549670733473680409194666521803651767 

1193123383965058588912264378749348068202790932902763246688455431*ci 

+4336608402971539347025237252135145455841035812345961936163213074515 

668913809930686458789021032092183097107068630380980654723975384806834 

399547265746928595076966467194028602153464501847919685450547334048541 

838603267464810745721708090718604776111130817877753185112995125825196 

9900571656166285598237774354646551889675112530575080545846831937*4 

-29667122481091026364901639098151640003965888696689225008930393424968 

602525252447427436681754274345497424043473832230915504701237012993339 

710812856574652111731271369111860472809751620103390836422947014561918 

596742170800703187841128631115509897386033559201960217645649903785096 

9895485234230652379590298859552730627182446420246826454883883689*4 

+8713668020397159975202209402993411331481527445489159653742816589285 

782666925647284363167235229251382288691898280884602035659684679625961 

618445548964513132925721174946043787600439882519732665839121162038859 

547938786254023880112832091100880202921766725381340769795219852772659 

31253523019517845556699435901635488455864886970958223801767254199*4 

-78305675171958244526488263512192141094477411091651712476670327217351 

782883871707407108769516955504305407699176179301277874612313054094467 
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597035177196760703985825236814433985312481352680953166404723733165962 
317279471444717834466650359718509004877505793494350037939473327528086 

0867972222715638415162988327352914888458003606019978160749679848*4 
-50522599203583621625595637115217448838607343066299095645637009117677 

845425111468649945267084093503724529301981623817384478696550110017337 
119404751862738670370910467517021683484411445974804570961028420156557 
202505669968683614281388897157201250167813030135004835015786338672614 
037510989474695713072179459238022542102052991592155908191171478*c| 
+35315707780022362747687233745930780256946702046940375643458709883815 
717091334531092687692414120601353572599659606631610854904141003914876 
375870980951760199224617620310327264953286419665073175139662877717976 
464011962351866981273603060052461996153251279426523988643704170224549 
3505041609179506395307679350757457119876008964575838569723696328*c9 
+26196064537923148987259844023868839991472305689605751242386858343654 
158387703055190547211883116255260309881480267447955878169277786134938 
237906112297170134620139058793181439502365320984021683968720825019581 
715365315571505452391632141389861154616280828170929049200565639622657 
40795249676124179031191642089155403460321065906554129943039359*C3 
The polynomials G'3,6(s6, C2) = 0,6*3,5(55,02) = 0, 63,1(51, C2) = have the 
same form as 62(02, C3) = given above. 

The polynomial 61(02) = has only four real roots, which are: —1.871987136, 
-0.3332612900, -0.2556946431, 0.2733107997; so the rank of the dataset is greater 
than 7. We embed it by joining the following vectors to the three shces: = 
(1 0), V2 = V3 = 0. The embedded dataset is = (X; vi)', X^ = (X'2 v^)' 
and X| = (X3 Vg)'. The rank of the embedded dataset will be calculated by two 
distinct methods. 

First, for the embedded dataset we see that ncolF =nrowr =8, so we can cal- 
culate the determinant of T as in Example 5, which is: 

det{T) =0 = 111296195967997*c^-163212875913821*ci-288078435761246*ci*C3 

+ 188384423078426*ci+139757151961919*cf C3-123835533958927*cfc^ 
+3188520736473*C2+1745777654358*C2*C3+145702375007129*C2*ci 
+ 154156258186696*C2*ci-30068441704134*C3-78231890782721*ci 
-9292669314727*c|+24148992371016*4 

Following the argument in Example 5, we note that there is a shght possibility 
that there will not be eight distinct values of (62, C3) such that the det{T) —0, because 
it is of degree 4. So, in general, following this approach of computing we can not assert 
that typical rank of generic 7x4x3 arrays is {7, 7 + 1} . However, let us continue 
our computation as in Example 5. We obtain the C matrix of Lemma 1, rounded to 
2 decimal digits. 
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[1,1,1,1,1,1,1,1] 
C I [5.68,63.84,0.6,44.31,3.49,9.33,3.93,21.29] 

[-38, 91, -1, 63, -23, -63, -26, 30] 
where the third row of C represents the randomly generated eight £3 values, and 
the second row represents the corresponding eight Ci values obtained by solving the 
(iet(r) =0 after plugging the £3 values in it. 

Now, we can obtain the S matrix of Lemma 1 by solving 



s'(X^ - c^X^) = 0' for A; = 2, 3 



:i7) 



eight times: The ith column of S corresponds to the eigenvector of the matrix 
r(c2i, c-zi) associated with the unique null eigenvalue 



S := 10-=^ 



X 



\ 



B := 10" 



X 



[-6.85, -6.54, -4.89, -6.51, 7.02, 6.74, 6.97, 6.44] 
[-4.02, 5.45, 8.94, 5.44, 3.85, 4.11, 3.90, 5.40] 
[9.29, -5.99, 14.8, -5.97, -9.49, -9.16, -9.43, -5.91] 
[-16.4, -6.87, -34.4, -6.83, 16.7, 16.2, 16.6, -6.71] 

[-4.31, 0.995, -6.31, 1.01, 4.15, 4.40, 4.20, 1.08] 
[-11.5, -5.37, -43.7, -5.36, 11.8, 11.4, 11.7, -5.32] 
[-3.22, -6.50, -6.24, -6.50, 3.20, 3.24, 3.20, -6.52] 
\ [-1000, -1000, -1000, -1000, -1000, -1000, -1000, -1000] / 
The matrix B' = S'Xi = A~^Xi of Lemma 1 is 

/ [1.06, -1.47, 22.7, -2.11, -1.84, -0.625, -1.61, -4.37] \ 
[-2.25, -1.24, -171, -1.79, 3.73, 1.35, 3.30, -3.72] 
[2.49, -0.077, -22.7, -0.111, -4.22, -1.48, -3.70, -0.228] 
\ [3.85, -0.280, -69.4, -4.00, -6.39, -2.31, -5.65, -0.806] / 

And A' = S"-*^; finally, we obtain X^ = AD(cj,)B' = X^. This was numerically 

verified. 

A second approach to compute the matrices S, A, B and C is via the Grobncr basis 
for the embedded system (17) using the lexicographic order (si, S2, S3, S4, S5, se, C3); 
note that ci is a free variable. The first Grobner basis polynomial is 

^1(02,03) = = 111296195967997*c|-163212875913821*c^-288078435761246*ci*C3 
+188384423078426*c2+139757151961919*c2*C3-123835533958927*cf ci 
+3188520736473*C2+1745777654358*C2*C3+145702375007129*C2*c2 
+ 154156258186696*C2*c|-30068441704134*C3-78231890782721*ci 
-9292669314727*c^+24148992371016*c^, 

which equals det(r). This shows that both approaches are identical for this par- 
ticular problem. 

4. Another Application of The Main Lemma 
Consider nongeneric dataset of size 4x4x3 
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/ [-872410,509152,-155756,301976] \ 
_ [-669515, 355308, -105576, 215236] 
^ [349983, -898362, 265770, -79182] 

\ [3285, -185950, 180998, 97398] / 
/ [-403995,481229,24054,201485] \ 
._ [-243133,337616,-4344,94484] 
^ [317091, -174294, -2454, -206076] 

\ [-317457,112640,183938,289254] / 
/ [-274447, 214327, -280750, 108851] \ 
_ [-252456,116912,-145020,92016] 
^ [-127464, -713802, 599526, 54318] 

\ [-38790,-204608,236662,21168] / 
To see if the rank of X is 4, we solve the system (9) composed of 8 polyno- 
mial equations in 5 variables via Grobner basis using the pure lexicograhic order 
(si, S2, S3, C3, C2). The elements of the Grobner basis are 

Gi{c2) = = -266104 + 1131869C2 + 1855673c^ - 10091484c3 + 3934656c^; 
^2(02,03) = = 70150154675210213-61657878275323159c2-700780737688415568( 
+30889123676791 1424c3 + 628616789525725c3; 
G'3,3(c2,S3) = = -79011958683266608845932181098557 

+37717083374737443006703954200886c2+901077269210427745705210304730192ci 
-390331538460948950190867958454016ci+1099565644871457602013702982455s3; 
^3,2(02,52) = = 39353103064214280234416428219949 

-l'90977009897456095062042069799807c2+209064381129999539517569775784236c 

-57940861139941085694575004742848ci + 5497828224357288010068514912275s2; 

G'3,i(c2,si) = = -29281575292540618957256320186316 

-3959967531501755611631756716147c2+390527303469098244882389504161956c2 
-165798278803217428934162052760128c^+1099565644871457602013702982455si. 
The polynomial Gi(c2) = is of degree 4 and it has four real roots, which are: 
-.3369565217, .2929292929, .2962962963, 2.312500000. So the rank of the dataset is 
4 by the main Lemma. Such datasets have been characterized by their defining 
equations in Landsberg and Manivel (2006). 

5. INDSCAL ARRAYS 

Let X e M^^^'^^^'^ be a tensor of order 3, where the ith slice Xj G M*^^*^ for i = 1, 
is symmetric. INDSCAL, proposed by Carroll and Chang (1970), is a statistical 
method used in psychometrics to analyse such arrays. For this reason, we shall 
name such an array an INDSCAL array to distinguish it from a general three-way 
array Y e ]gJxJxK discussed above, where such a decomposition is usually named 
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PARAFAC, see Harshman (1970). A rank 1 INDSCAL array or a decomposed tensor 
is 

D = a«)b«)b, (18) 

where a G and b G M"^. 

The following theoretical results are known for generic INDSCAL data X e 
^/xjxj . gy Zellini (1979), sec also Rocci and Ten Berge (1994), if 7 > J(J + 
l)/2, then rank(X) = J(J + l)/2. b) / x 2 x 2 and 1x3x3 arrays are studied by 
Ten Berge, Sidiropoulos and Rocci (2004). The rank computation problem has also 
been approached from a numerical point of view by Comon and ten Berge (2008), 
who applied applied Terracini's lemma, based on the numerical calculation of the 
maximal rank of the Jacobian matrix of (2), to obtain numerically the generic rank 
of some INDSCAL three-way arrays. The numerical method based on Terracini's 
lemma, when used to evaluate rank over R, gives the generic rank when the typical 
rank is single- valued, and the smallest typical rank value otherwise. 

For INDSCAL data (7) becomes 

s'Xfe = bkh' for k = 1,..., J, (19) 

for Xfc e W'^, s e and b e 

We note that (19) has two indeterminacies, scale and sign: It can be rewritten as 

s^'Xfc = bkh' for k = 1, J, , where for instance, = As for any scalar A > and 
b = A^/^b. The second indeterminacy is the sign indeterminacy of b : replacing b by 
— b in (19) does not change the equality in (19). To eliminate both indeterminacies, 
hereafter, we fix 

hi = 1. (20) 

We will represent the set of solutions of (19) subject to (20) by V (Veronese 
variety) . 

We are interested in the study of the number of solutions of (19) subject to 
(20) over R for genenc INDSCAL data for 2 < J, / < J(J + l)/2. We distinguish 
three cases named, minimal when / = l-(-J(J — l)/2, overdetermined when / > 
1-|-J(J — l)/2, and, underdetermined when I < 1 + J{J — l)/2. The overdetermined 
systems is similar to the one discussed above. 

Theorem 3 (minimal system= Veronese variety): Let 7=1 + J(J — l)/2 
and 2 < J < 7, then for generic INDSCAL data the number of roots (real or complex) 
of the polynomial system (19) is 



(21) 
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Proof: Let [bi, be an element of the projective space P'^'^^^K We note that the 

right hand side of (19) is a Veronese variety of degree d = 2, which is the image of 
the Veronese map, 1/2, defined by 



by sending 



[61, bj] [bl, 6162, Mj-1, 

where the image has N + 1 — elements composed of binomials in bi,...,bj. 

While the left hand side of (19) is a general linear space of projective dimension 7 — 1. 
The number of intersections of the general linear space with the Veronese varity is 
finite, when / — 1 = A^— (J — 1); that is 

7 = 1 + J(J- l)/2. (22) 

When (22) is true, the finite number of intersections is the degree of the Veronese 
variety given in (21), see for instance Harris (1992, p. 231). 

Corollary 1: The typical rank of INDSCAL arrays with a minimal system have 
more than one rank value and the minimum attained value is 7. 

Proof: For minimal systems and 2 < J < I, I < degV. The rank of a generic 
INDSCAL array with a minimal system is 7, if the number of real roots of V is 
greater than or equal to 7; otherwise its rank is greater than 7. 

5.1. Example 7. We consider a simulated generic dataset of size 4x3x3 having 
the following 4 slices 

[54, 107, 161] \ / [114, -49, -125] 

Xi := I [107,58,13] X2 := [-49,-144,-76] 
[161,13,134] / \ [-125,-76,-8] 

[-44,7,-48] \ / [50,92,-4] 

X3 := I [7,-36.-11] X4 := [92.100,1] 

-48,-11,-154] / \ [-4,1,-100]] 

INDSCAL 4x3x3 arrays have been studied in detail by Ten Berge, Sidiropoulos, 
and Rocci (2004), where it is shown that if a certain polynomial of degree 4 has 4 
real roots, then rank(X.) — 4, otherwise the rank is 5. 

The Grobner basis with pure lexicographic order given by the following sequence 
(61, 62, Si, S2, S3, S4) of the free variables is formed of 6 polynomials listed below. The 
first polynomial 6*4(54) = is of degree 4, as shown by ten Berge, Sidiropoulos and 
Rocci (2004) and Theorem 3, and it has 2 real roots -0.1881015674e-2, 0.7632125093e- 
1, so the rank of the dataset is greater than 4. 
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6-4(54) =7337669360341773654444527-4293727819369270858661345768*54 
-211863296775796994233864209576*s|-1486920579131214046506874714272*s3 

+65728692033647334980166673748496*5^ 

^3(53,54) =17560973913802573904674715803175900627366683113948054931 

-161223257377160952551452668901669378891965735728005400638*54 

-5837696072380594108159410240186154115595431615363079926796*5^ 

+1318397923701745931624444235931465979525756973974101183796472*5| 

+5238806078525191567165441234720094289579153419952152373008*53, 
^2(52,54) =9628239825303370993360207993471191965769478430007385299 

+11068252823433558277754489464864186098568340630108319931766*54 

+705995008022931051200292932727627624011448550604188031890116*5^ 

-9118208129962992736609222153263187238797217283937992554624760*5| 

+20955224314100766268661764938880377158316613679808609492032*52, 

Ci (si , 54) =-9384923940854434492010502183000235854601288296806481877 

-1964822700901995622589398714750515451903512627004147123640*54 

+1914122466041979887104614509203684754556718777770379683036*5^ 

+689540537276652567152783462787012635706612460537819392520176*5^ 

+2619403039262595783582720617360047144789576709976076186504*51, 

G'2,4(&2,S4) =-155028823048701914654384720617223421617571775035134363431 

-66211886029016478638439782073183667186014332650277499221046*54 

-3324980827880602990137773057985895339331874506372137213739628*5^ 

+44167907819442655873419676087304190072683931413772512430408456*5^ 

+1309701519631297891791360308680023572394788354988038093252*^2, 

61(61,54) =-1514819909584866108143372567736179608809277026174154396973 

-345404316274377016240902795956549557726567542396897273802586*54 

-7696649874609748839698115121543942313397062924509594873506300*52 

+161342893355178852699731325084429928324620333696539365418650824*5^ 

+1905020392190978751696524085352761559846964879982600862912*6i 

5.2. Example 8. We consider a simulated generic INDSCAL dataset of size 7 x 
4x4 having the following 7 slices 



Xi : = 



X, := 



/ [140,86,-110,-4] \ 
[86, -182, 70, 36] 
-110,70,104,183] 
\ [-4,36,183,148] / 
/ [178,15,-186,52] \ 
[15,196,119,-148] 
[-186,119,-138,43] 
\ [52,-148,43,-110] / 



X, 



X4 := 



/ [-20, 100, 173, -56] \ 
[100, 128, 101, 75] 
[173,101,124,65] 
V [-56,75,65,-158] / 
/ [-8,-137,21,20] \ 
[-137,-60,64,5] 
[21,64,-24,-14] 
V [20,5,-14,-128] ) 
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/ [194,-164,-36,-6] \ 



/ [-22,74,85,-40] \ 
[74, -198, 23, -53] 
[85,23,-152,18] 
\ [-40, -53, 18, -48] / 



[-164,2,-114,-110] 



[-36,-114,64,-127] 

V [-6,-110,-127,-18] / 
/ [-94,109,-16,90] \ 



[109, 124, 164, -93] 
[-16,164,98,-134] 
\ [90,-93,-134,184] / 



The Grobner basis with pure lexicographic order given by the following sequence 
(si, S7, 6i, 62, ^3) of the free variables is formed of 10 polynomials, but only the 
first one is shown below. The first polynomial ^3(63) = is of degree 8, as shown 
in Theorem 3, and it has 2 real roots -4.615952848, 1.035693119, so the rank of the 
dataset is greater than 7. 



Gsih) =-267319790697212354162205439965563724346086890209668628287611296 

-418573483979735109514695930195818332286961955303805928337210144*63 

-53224562968122644847846329140305933491773608156555442814188832*6| 

+190260260230311283947025128614232688395842607775283676963072480*6^ 

+172806709139583658797792038234309205181892588602072553465939944*61 

+164461253658569745584828839860332360925297521683057843681458288*61 

+95090716874891491062104108972298040778579595301835996945880112*6^ 

+24575461542028106507015735163996805821952192308876843008456228*61 

+2240887382441309183839416634048576470976843441637962999441259*61 



We introduced a new method to compute ranks of three-way arrays, by showing that 
it is intimately related to the solution set of a system of polynomial equations, which 
is a well developed and active area of mathematics known as algebraic geometry. The 
new method was used to compute numerically the ranks of some sizes of three-way 
arrays over R via Grobner basis. 

The problem of computing the rank of overdetermined systems by solving embed- 
ded polynomial systems is a work in progress. 
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Appendix 1 

Below the matrix YA; = X'^.. 
> K := 3; 



6. Conclusion 
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> J := 8; 

> I := 15; 

> with(LinearAlgebra); 

> Yl := RandomMatrix(J, I); 

> Y2 := RandomMatrix(J, I); 

> Y3 := RandomMatrix(J, I); 

> S := Vector(l .. I, 1); 

> for h from 1 to I-l do 
S[h] := Sh end do; 

> Ml := Y2-C2*Y1; 

> PI := Ml.S; 

> M2 := Y3-C3*Y1; 

> P2 := M2.S; 

> poly := [seq(Pl[l], 1 = 1.. J), seq(P2[n], n = 1 .. J)]; 

> witli(Groebner); 

> liste := seq(s/i, h = 1 .. I-l); 

> polyG := Basis(poly, plex(liste, C3, C2)); 

> polyG[l]; 

> S := [solve(polyG[l])]; 

> nops(S); 

> fsolve(polyG[l]); 

> nops([fsolve(polyG[l])]); 

> fsolve(polyG[l], a, complex); 

> nops([fsolve(polyG[l], C3, complex)]); 
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